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ABSTRACT 

We investigate the use of Antithetic Variables, Control Vari- 
ates and Importance Sampling to reduce the statistical errors 
of option sensitivities calculated with the Likelihood Ratio 
Method in Monte Carlo. We show how Antithetic Vari- 
ables solve the well-known problem of the divergence of 
the variance of Delta for short maturities and small volatil- 
ities. With numerical examples within a Gaussian Copula 
framework, we show how simple Control Variates and Im- 
portance Sampling strategies provide computational savings 
up to several orders of magnitude. 

1 INTRODUCTION 

Monte Carlo (MC) simulations are one the main tools em- 
ployed in the Financial Services industry for pricing and 
hedging derivatives securities. In fact, as a result of the 
ever increasing level of sophistication of the financial mar- 
kets, a considerable fraction of the pricing models employed 
by investment firms is too complex to be treated by ana- 
lytic or deterministic numerical methods. For these models, 
MC simulation is the only computationally feasible pricing 
method. 

The main drawback of MC methods is that they are gen- 
erally computationally expensive. These efficiency issues 
become even more dramatic when MC simulations are used 
for the calculation of price sensitivities, i.e., the derivatives 
of the option price with respect to the parameters of the un- 
derlying model, also known as Greeks. In fact, the standard 
method for the calculation of a price sensitivity, say with 
respect to a parameter 6^, is based on a finite difference 
approximation of the derivatives dV(9)/d9i ( . This method, 
also known as 'bump and reval', involves repeating the MC 
simulation, and evaluating the finite difference estimate 



dv(o) v(e u ...,o k +AO,...,e n )-v(e) 



do, 



AO 



(1) 



for a small increment AO. The main virtue of this method 
is that it is straightforward to understand, and it requires 
minimal implementation effort. The drawback is that ad- 
ditional MC simulations are required for each sensitivity, 
and that the finite differences (Q]) may be affected by large 
statistical errors, especially for payout with discontinuities 
(Glasserman 2004). As a result, hedging derivative securi- 
ties with MC simulations can be extremely time consuming. 

Alternative methods for the calculation of price sensi- 
tivities have been proposed in the literature (for a review see 
Glasserm an 20041 1. Here we concentrate on the so-called 
Likelihood Ratio Method (LRM). The principal advantage 
of this technnique when compared to 'bump and reval' is that 
it allows to calculate all the sensitivities simultaneously in a 
single MC simulation, and a single set of payout evaluations. 
In addition, the variance properties of LRM estimators are 
not affected as much by discontinuities in the payoff. As a 
result, for digital and barrier options, LRM may provide a 
better convergence than bumping (Glasserman 2004). The 
main drawback is that the statistical uncertainties of LRM 
estimators are nonetheless difficult to predict, and can be 
sometimes large. What is worse, in some cases such un- 
certainties are even known to diverge thus making the MC 
simulation very time consuming if not hopeless in practice. 

In order to address this difficulty, in this paper we 
investigate three Variance Reduction techniques - Antithetic 
Variables, Control Variates, and Importance Sampling - that 
can dramatically improve the MC convergence of LRM 
estimators. In the next Section, we begin by reviewing the 
rationale of LRM, specializing our discussion to a Gaussian 
Copula framework very common in the financial practice. 
The use of Antithetic Variables is discussed in Section[3] In 
particular, we will show how this simple technique solves 
the well-known problem of the divergence of the variance of 
LRM Deltas for short maturity and low volatility. Then, in 
Section|4]we illustrate how Control Variates and Importance 
Sampling can drastically suppress the statistical uncertainties 
of the LRM estimators thus reducing the computational cost 
for the Greeks by orders of magnitude. 
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2 LIKELIHOOD RATIO METHOD 

The arbitrage-free price of a derivative security can be 
expressed in general as the expectation value of the dis- 
counted cash flows, G(x), over a risk-neutral probability 



density (Harrison and Krep s 1979| , Pg{x), 



V(d)=Ep[G(x)} = / dxG{x)P e {x) 



(2) 



where x=(x\,... ,Xti) is a A^-dimensional vector represent- 
ing the underlying random factors upon which the claim 
is contingent. Here the vector 9 = (9i,...,9 n ) represents 
a set of parameters whose value is generally determined 
by calibrating the chosen model or, equivalently, the den- 
sity Pq(x), on the prices of securities liquidly traded in the 
market. 

Whenever the dimension N of the state variable x is 
large (say N >4) MC methods are the only feasible route 
for estimating expectation values of the form (0. In their 
simplest incarnation, these consist in averaging the payout 
function G(x) over N p independent random realizations of 
the vector x, say x[m], generated according to the probability 
density Pg(x), 
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Figure 1: Delta of a European Call option (f2Tb for So — 100, 
K = 25, and r = 0.05 as a function of the time to maturity 
T for a = 0.1 (a), and as a function of the volatility a for 
T = 1/12 (b): crude MC (triangles and continuous line), 
MC with Antithetic Variables (circles and dashed line). 



V{9)~V 



— £ G(x[m] 



N, 



(3) 



In fact, the central limit theorem ( Kallenber g 1997[ ) ensures 
that, for big enough samples, the values of the estimator V 
are normally distributed around the true value, and converge 
for N p — > °° towards V namely 



V^-^G(x[m})±- 



INn 



(4) 



where fc 2 = Ep [G(jc) 2 ] — Ep [G(jc)] 2 is the variance of the 
MC estimator. Here, following a common terminology, 
we refer to K"/ A 



/Np as statistical uncertainty or statistical 
error. Although Eq. (0J ensures the convergence of the 
MC average to the expectation value (ffj) provided that K is 
finite, the square root law in can make the calculation 
of accurate estimates time consuming. 

This is particularly true, for the MC calculation of the 
Greeks. In fact, the variance and bias properties of finite 
difference estimators of the form (fl} can be in some cases 
rather poor (Glasserman 2004). This is because, while the 
bias of the finite difference ([TJ can be made in general arbi- 
trarily small by reducing the value of A9, its statistical error 
can in some common cases diverge for A9 — > 0. When this 
happens, choosing a value of A9 small enough to reduce the 
bias to an acceptable level may require a large computational 



cost in order to obtain statistically accurate results. For a con- 



vergence analysis see e.g., (Milstein and Tretyakov 2005). 

Several methods have been recently proposed in the 
literature in order to speed up the calculation of option 
sensitivities (Glasserman 20041) . Here we will concentrate 
on the so-called Likelihood Ratio Method (LRM). Under 
mild regularity conditions on the probability density Pg(x), 
the sensitivity of the option price (f2]i with respect to any 
parameter 9k can be obtained as 



Qk = mp-=MG(x)n k {x) 

o9 t 



(5) 



i.e., by calculating the expectation value of the original 
payout function multiplied by the so-called Likelihood Ratio 
weight 



giving as MC estimator: 
1 



d\ogP {x) 
d9 k 



Qk- Jf £ G(x[m])Q k (x[m]) 



(6) 



(7) 



In the following, we will specialize our discussion to 
the case where the probability distribution Pg(x) is a N- 
dimensional Gaussian Copula. The latter is defined by a 
correlation matrix E, and a set of marginal cumulative 
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distributions M,-(xj), so that the joint distribution reads 



F(x) = f[ /' d yi P(y u ...,y N ) 

= <t> N (4>- 1 (Mi (*, )),.,*- 1 (M w (x w )) ; E) . (8) 

Here $/y(jCi,... is the cumulative distribution of a 

Af-dimensional Gaussian random variable with zero mean 
and correlation matrix E, and <t>(.) is the standard normal 
cumulative distribution. 

If we restrict to the case in which the correlation matrix 
E does not depend on the parameters 9, the LRM weight © 
for the Gaussian Copula distribution © can be expressed 
as 

N 

Cl e (x) = £a e logm,(x I )-Z(x) r (E- 1 -I)d B Z(x) , (9) 

i=i 

where / is the A^-dimensional identity matrix, m,(^/) = 
dM{xi) I 'dx[ is the r'-th marginal density function, and the 
components of the vector Z(x) are 



Z, =*- 1 (M,-(x i -)), 



and those of dgZ(x) read 



d e Zi = 



deMi(xi) 



(10) 



(11) 



with (j>(.) the standard normal density function. A derivation 
of Eq. (0 is given in the Appendix. 

Gaussian Copula models of the form (0 are widely 
used in Financial Engineering. Indeed, the above formu- 
lation can be used to evaluate structured European options 
written on several assets, e.g., equity, commodity, rates 
or foreign exchange pairs (IHull 20021 1 . In this case, the 
marginal distribution of each asset is typically implied from 
liquidly traded Vanilla options, and the co-dependence be- 
tween the factors is modeled by means of the Gaussian 
Copula. Structured credit pricing, e.g., for CDO and CDO 2 
dSchonbucher 20 03). can be also performed within a similar 
framework. In general, whenever the marginal distributions 
above are not known in closed form, e.g., they are calculated 
numerically by means of a calibration procedure, the deriva- 
tives in Eq. (fTTT i can be easily computed by means of finite 
differences. This does not generally introduce accuracy 
or stability problems provided the calibration algorithms 
employed are numerically stable. 

It is easy to see that the LRM weights above give the 
expected result in the case a multi-asset lognormal model 
of the form 

S i = $ex V [( r -o?/2)T + a i VTZ i ] , (12) 



where Sf and o, are the spot price and volatility of the ;-th 
asset, T is the maturity of the option, and Z, are standard 
normal variables with correlation £y = E[Z,-Z ; -]. In fact, one 
clearly has Xi = Zi with 



Z 



logSi/$-{r-o?/2)T 

OfVr 



(13) 



so that the LRM weights for the z'-th Delta and Vega 
dHull 20021 1 read respectively 



E-iZ]; 



mz)=^o g p ei x)=^ 



(14) 



and 



Q.Uz) = ^-logP (x) = (^--Vf)[Z- l Z] i -- . (15) 

dOi \C7, / Gi- 

lt is straightforward to verify that these equations are in 
agreement with the general expression for the LRM weight 
in a Gaussian model given in (Glasserman 2004) 



a e (Z) 



■~Tr[tr 1 det]+~xZ- i (d e t)t- 1 X 



X±~ l d G m 



(16) 



with nit = log 5? + (r - a 2 /2)T, X t = (JiVTZi, and E, 7 = 
CT;CT/E, 7 . 

In the special case of a single asset, the weights above 
simplify to the well-known expressions 



n A (z) = 



aVfS ' 



and 



Qv(Z) 



Z l -\ 



-zVf 



(17) 



(18) 



3 SOLVING THE PROBLEM OF DELTA'S 
DIVERGING VARIANCE WITH ANTITHETIC 
VARIABLES 

For a given number of MC iterations the calculation of the 
Greeks by means of LRM is generally fast when compared 
to bumping. However, the speed of convergence of the 
LRM estimators is difficult to predict a priori for a given 
problem, it is payout and parameters dependent, and can 
be in some instances particularly slow (Glasserman 2004). 
In fact, since the LRM weight has in general zero mean as 
a result of the identity 



d e J dxP e (x)=0 



(19) 
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the LRM estimators for the Greeks have no definite sign. 
This can give rise to poor variance properties whenever the 
configurations with opposite sign have similar weight in the 
MC average © so that the final outcome is the result of the 
cancellation of two comparable and not necessarily highly 
correlated quantities. 

Table 1: Delta of the Basket Call option d27t f or T = 0.5 
and r = 0.05. The volatilities of the assets are all equal to 
o = 0.3, and their Forwards range between 51.3 and 55.9. 
The uncertainties are reported in parenthesis. 
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Figure 2: Vega of the Basket Call option (|27T i as a func- 
tion of the strike price. Crude MC: triangles and dashed 
line. MC with Variance Reduction techniques: circles and 
continuous line. Empty circles: Antithetic Variables with 
Control Variates. Full Circles: Antithetic Variables with 
LSIS. 



In particular, a common problem generally reported 
in the literature (Glasserma n 2004] IJaeckel 20021 ) is the 
divergence of the variance of the LRM weight for Delta 
Eqs. ( TT4l and ( TTTb in the limit of small volatility and short 
maturity. Indeed, from Eq. dT4b 

Var[£2y oc J_ - oo (20) 

for (JiVT — > 0. This is illustrated in Figure 1 for a 'deep in 
the money' Call option (i.e., with very low strike compared 
to the expectation value of the underlying asset, or Forward, 
see lHull 20 02) on a single lognormal asset with undiscounted 
payout 

G(S) = (JS-K)+ (21) 

where K is the strike price. The LRM estimator becomes 
extremely noisy for T — > and for small volatility, to the 
point of not providing any useful information for maturities 
shorter than a few weeks for any practical number of MC 
iterations. 

The divergence of the LRM weight for Delta is due to the 
break down of the absolute continuity property of the proba- 
bility density function, which is required to take the deriva- 
tive inside the expectation in Eq. © (Glasserman 2004). 
This generally affects the LRM estimators (O for Delta also 
for non-lognormal models. 

Although, to the best of our knowledge, it was not 
previously noted in the literature, this problem can be easily 
overcome by using Antithetic Variables. Indeed, since the 
LRM weight (fT4l is odd in each of the Gaussian random 
variable Z,, the Antithetic estimator (Glasserman 20041 for 



the weight reads 

ftkU = M ' 2 AV = (22) 

so that Ep[Q^|„ nf ] =0, with zero variance. What is more, 
it is also possible show that the variance of the product of 
the payout and the weight in ^ is generally bounded as 
(TV? -> 0. This can be realized by means of the following 
simple heuristic argument. Consider for simplicity a single 
asset payout under the lognormal model ( fl~2~b . This can be 
approximated for small maturities and volatilities as 

P{Z) ~c + Cl aVfz + O{a 2 T) (23) 

with co and c\ constants . As a result, the LRM Antithetic 
estimator for Delta (0 reads in this limit 

QaU - c 2 ^ s0 + ci + 0(aVf) (24) 

whose variance is clearly bounded for 0\/T — > as a result 
of the cancellation of the leading term. 

Figure [T] illustrates the efficacy of the method: the 
Antithetic LRM estimator provides a stable Delta and a 
practically constant statistical uncertainty for all maturities 
and volatilities. The variance reduction with respect to the 
crude LRM estimators is around 10 for 1 year, and around 
500 for 1 week maturity. 

While the results presented here are for a lognormal 
model, we have found comparable variance reductions also 
for market-implied marginal distributions M, (x) for a variety 
of equity, foreign exchange, interest rate and commodity 
underlying assets. Indeed, the antisymmetry of the LRM 



Luca Capriotti 



weight (O for Delta is generally satisfied, at least approx- 
imately, also for skewed distributions. It is also worth 
noting that although here we have limited our discussions 
to European-style payouts in a Gaussian Copula framework, 
it is easy to realize that Antithetic Variables generally solve 
the problem of the divergence of the variance of Delta esti- 
mators also for path-dependent options (similar ideas have 
been also used by Mike Giles in the implementation of the 
so-called Vibrato Monte Carlo technique (IGiles 20071 0. 

4 CONTROL VARIATES AND IMPORTANCE 
SAMPLING 

Although the variance of LRM estimators generally remains 
finite when using Antithetic Variables, their convergence 
can still be poor in some cases, especially for Vega. For 
this reason, an efficient implementation of LRM generally 
requires the use of other Variance Reduction strategies. 
Here we consider Control Variates and Importance Sampling 
dGlasserman 2004l >. 

The Control Variates method aims to reduce the statis- 
tical uncertainty of a MC average by exploiting the corre- 
lation of its statistical samples with those of some quantity 
(the control) whose expectation value is known a priori 
(Glasserman 2004). This technique may result in spectac- 
ular variance reductions but requires some closely related 
estimator with a known integral over the sampled probability 
distribution. 

For this discussion, we will restrict ourselves to controls 
that in our context are generally readily available, and we 
will consider the LRM weights and the derivatives of the first 
moment of the marginal distributions in Eq. ||8). In fact, as 
previously mentioned, the LRM weight © has always zero 
expectation value. It is also reasonable to expect the weight 
to be somewhat correlated with the corresponding Greek 
estimator in (0. On the other hand, the first moment of 
each marginal distribution is usually known as it contains 
the information on the Forward of the underlying asset, 
which is the first thing to be usually matched with the 
available market data. Its derivatives with respect to the 
main model parameters are also generally known. For 
instance, in a typical financial context, including models 
with skew dHull 20021) . one has 



dE[Sj~ 

~Jsf 



P(0,T) 



(25) 



where P(Q, T) is today's (t — 0) price of a zero coupon bond 
maturing at time T, and 



dE[Sj] 
dcTi 



= 



(26) 



In the following, when calculating the Greeks of the i-th 
asset, we will use the corresponding LRM weight, and Risk 
of the Forward. Note, however, that when calculating the 
LRM Delta with Antithetic Variables, we can only use the 
Delta of the Forward as a control, as the Antithetic estimator 
for the LRM weight (f22j) is identically zero. 

Table 2: Same of Table Q] for Vega. 
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Importance Sampling techniques, on the other hand, do 
not rely on the knowledge of any closely correlated estimator 
but aim to reduce the variance by sampling more effectively 
the domain of integration in Eq. (f2|l. Here we will use a 
recently introduced Importance Sampling strategy based on 
a Least-squares optimization, namely the Least-Squares Im- 



portance Sampling (LSIS) ( [Capriotti 2008| |Capriotti 2 007). 
In particular we will use as trial densities mean-shifted sin- 
gle mode and bi-mode multivariate Gaussian distributions 
(see references above). 

In Table [T] we compare the results obtained by using 
Antithetic Variables only, and in combination with Control 
Variates and Importance Sampling for the Delta of one of 
A' = 10 assets of a Basket Call option with undiscounted 
payout 



G(S) 



N 

(=i 



K 



(27) 



in a lognormal model of the form dT2l ). Here, as an indicator 
of the efficiency gains introduced by the different methods, 
we have defined the variance (efficiency) ratio as 



ER 



a[CrudeMC] 
ofVR] 



(28) 



where the numerator and denominator are respectively the 
statistical errors (for the same number of MC paths) of 
the crude estimator and of the one obtained with different 
Variance Reduction techniques. 

As shown in Table [T] for the considered maturity, the 
efficacy of Antithetic Variables decreases moving away from 
the deep in the money region, and is generally limited to 
around a factor of two for larger strikes. Using the Delta of 
the corresponding Forward $25[ as Control Variate provides 
about an order of magnitude variance reduction around the 
'at the money' point (i.e., for strikes around the Forward 
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value). Decreasing further the strike results in a spectacular 
suppression of the statistical errors of the Control Variate 
estimator. This is expected due to the very high correlation 
between the option payoff in this regime and a simple 
Forward contract. However, as the strike moves instead in 
the 'out of the money' region (high strikes compared to the 
Forward) such correlation rapidly decreases, and Control 
Variates become practically ineffective. 

On the other hand, a simple Importance Sampling strat- 
egy based on a single mean-shifted Gaussian trial density 
and implemented by means of the LSIS approach proves 
to be very effective at all level of moneyness. In fact, al- 
though Importance Sampling gives smaller efficiency gains 
than Control Variates for deep in the money options, it 
provides at least two orders of magnitude speed up across 
all strikes, including the out of the money regions where 
Control Variates lose their efficacy. Indeed, as it can be 
generally expected when using mean-shifted trial densi- 



ties (Capriotti 2008), Importance Sampling is particularly 
effective moving away from the money in either direction. 

LRM estimators for Vega are generally noisier than the 
ones for Delta. This can be understood from the form of 
the weight in Eq. (fl~5T > as it involves the second moments 
of the Random increments which have generally a larger 
variance. This is illustrated in Fig. [2] where we plot one of 
the Vegas for the Basket option above as a function of the 
strike: for low strikes, as Vega becomes smaller, the crude 
MC estimate becomes extremely noisy. As mentioned, this 
is due to fact that the small value of Vega in this regime 
is the result of the cancellation of two poorly correlated 
stochastic quantities representing the averages of the LRM 
estimator over the configurations in which it has a definite 
(positive and negative) sign. 

As expected, due to the presence in the LRM weight 
of even terms in the random increments, the efficacy of 
Antithetic Variables is very limited for Vega (see Table [2j. 
In contrast, using both the LRM weight and the Vega of 
the Forward as Control Variates provides sizable variance 
reductions both for in the money and at the money options, 
with the expected remarkable efficiency gains for small 
strikes. On the other hand, Importance Sampling - although 
producing smaller variance reductions for in the money 
options - becomes particularly effective for larger strikes. 

In this sense, Importance Sampling and Control Variates 
appear to be somewhat complementary for this problem, and 
Figure [2] displays the LRM Vega as obtained by choosing 
the most effective of the two techniques for each strike. 
This results in orders of magnitude savings in computer 
time with respect to the crude MC calculation. 

ACKNOWLEDGMENTS 

It is a pleasure to acknowledge useful discussions with 
Mike Giles, Paul Glasserman, David Shorthouse, Jacky 



Lee, Anton Merlushkin, Mark Stedman, and Sanjay Chawla. 
The opinions and views expressed in this paper are uniquely 
those of the author, and do not necessarily represent those 
of Credit Suisse Group. 

A APPENDIX: DERIVATION OF EQ. ® 

First, by differentiating the joint cumulative distribution <[8J 
one obtains the corresponding probability density function 



P(x) = <j> N (<t>- 1 (Mi(x l )),...,<S>- 1 (M N ( XN ))-X) 



x n 



^-\Mi{ Xi ))) 



(29) 



where m,(x,) ~dM(xi)jdxi is the 2-th marginal density 
function, and $n(xi, . . . ,xyv;E) is the multivariate Gaussian 
density with correlation E. Then, taking the logarithm of 
Eq. d29l gives 



logP(x) = £ (logm,-(*/) -log0(Z,-)) - l -Z{x) J 'EZ(x) 

i= 1 " 



ylog27T-ilog(detE) 



(30) 



where we have used the explicit form of 0at(Zi, . . .Z^;E), 
and the definition in Eq. (TTOb . Hence, the derivative with 
respect to 9 of the latter equation, when <9gE = 0, can be 
written as in Eq. (0 where dgZ(x) (TTTT > can be obtained by 
deriving Eq. (TTOb . 
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